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We study the transient heat current out of a confined electron system into a weakly coupled 
electrode in response to a voltage switch. We show that the decay of the Coulomb interaction 
energy for this repulsive system exhibits signatures of electron-electron attraction, and is governed 
by an interaction-independent rate. This can only be understood from a general duality that relates 
the non-unitary evolution of a quantum system to that of a dual model with inverted energies. 
Deriving from the fermion-parity superselection postulate, this duality applies to a large class of 
open systems. 
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Energy decay due to heat currents is of key importance 
in the continued downscaling of electronic devices [1]. 
The quantum [2-6] and interaction effects [7-9] that 
arise on the nanoscale give rise to new possibilities [10- 
13] and motivate both fundamental [14] and applica¬ 
tion oriented [3, 4, 8, 13, 15] studies on quantum heat- 
engines, possibly realized in, e.g., cold atoms, trapped 
ions, or quantum dots. The successful control and ex¬ 
ploitation of heat in nanodevices requires both a funda¬ 
mental understanding and the practical ability to detect 
and manipulate few-electron heat currents. Under sta¬ 
tionary conditions, progress has been achieved using var¬ 
ious approaches [16-18], including heat transfer through 
molecular-scale devices [19] with electrostatic gating [20]. 
However, any device is eventually adjusted by some ex¬ 
ternal agent that provokes a time-dependent response. 
In the context of electronic heat currents, this raises a 
very basic question that, despite recent promising theo¬ 
retical [21-27] and experimental [28-30] studies, has not 
been answered so far: how does a small electron sys¬ 
tem, typically governed by a strong level-quantization 
and Coulomb interaction, dissipate in time its stored en¬ 
ergy into a coupled electronic bath? 

The essence of time-dependent transport in such sys¬ 
tems is already captured by the simple model sketched 
in Figs. l(a,b). Here, an instant energy shift of a sin¬ 
gle electronic orbital in a quantum dot leads to a time- 
dependent charge current / 7 v(t) [31, 32] and heat current 
Iqit) into a tunnel-coupled electrode. In the weak cou¬ 
pling regime, expressions for these currents can be calcu¬ 
lated straightforwardly, and in the case of the transient 
charge current I Nit) also allow for an intuitive physical 
understanding [33]. This is, however, not the case for the 
heat current Igit) = ace“'''‘=‘4-0^6“'’'*’*. Compared to the 
charge current I nit) oc , the heat current contains 
a second decay mode. The mere presence of this mode 
can be expected: it originates from the dissipation of the 
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FIG. 1. (a,b): Quantum dot with charging energy U sub¬ 

ject to a gate-voltage switch. The instant level shift from 
£0 to e causes two electrons to be sequentially expelled, (c): 
Inverted stationary state for the dual model obtained from 
the final state (b) by inverting all energies, as prescribed 
by Eq. (4). (d,e); Amplitude of the charge mode, Uc and 
of the parity mode, Op in the time-dependent heat current 
Igit) = OcC -I- OpC 'r*, plotted versus the initial (eo) and 
final (e) level position for T = 0.11/ S> F. The black crosses 
mark the switch in (a,b), the circles a switch to e = fi. 


Coulomb energy. However, what is quite remarkable is 
that its rate 7 p turns out to be completely independent 
of the interaction strength U [34-37] - despite entering 
the heat current only as a consequence of the interaction. 
Even more surprisingly, as indicated by the blue dashed 
line in Fig. 1(e), both excitation amplitudes, Oc and Up, 
show an abrupt change at an energy typically associated 
with electron-electron attraction [38] , even though we are 
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dealing with a system governed by repulsive interaction. 
As we will show, both surprising effects have a fundamen¬ 
tal origin and they can occur in a broad class of systems. 

The difficulty in explaining the counter-intuitive ob¬ 
servations in the heat current ultimately stems from the 
non-unitarity of the dynamics of open quantum systems. 
In closed systems, the hermiticity of the Hamiltonian im¬ 
plies that its left and right eigenvectors are simply ad¬ 
joint to each other. Thanks to this property, the ampli¬ 
tude with which an energy eigenstate enters the full time 
evolution of the system is simply given by the overlap 
between this eigenstate and the initially prepared state. 
For dissipative open systems, a similar and equally in¬ 
sightful relation for the amplitndes of the decay modes 
is not known. Yet, such a relation would be extremely 
valuable, not only for the understanding of the tran¬ 
sient heat-current in Fig. 1, but also for a broad class of 
other dynamical problems, ranging from qubit dynamics 
to molecular electron-transfer. 

In this Rapid Communication, we identify a duality 
relation between decay modes and amplitudes. This du¬ 
ality applies to all open fermionic systems that can be 
modeled by a Hamiltonian = H + -|- with an 

arbitrary local Hamiltonian H and (a) effectively non¬ 
interacting, wide-band fermionic reservoirs H^, with a 
coupling that is (b) bilinear in the fermion fields 
and (c) energy-independent, but otherwise arbitrarily 
strong. Under these very general assumptions, the dual¬ 
ity emerges mainly as a consequence of the fermion-parity 
superselection postulate of quantum mechanics [39-42], 
by which states with even and odd fermion number can¬ 
not be superposed. Applied to the system shown in 
Fig. 1, it explains in full detail the poorly understood 
heat decay. 

Mode-amplitude duality for open systems. An open 
quantum system is most naturally described in terms of 
its reduced density matrix p(t), whose equation of mo¬ 
tion is dtp{t) = —i[H,p{t)] + J*dt'yV{t — t')p{t') [with 
h = kB = e = 1]. The kernel W takes into account the 
coupling to the external reservoirs [43-46] that causes 
p{t) to decay. Introducing the Laplace transform of 
yy, W{ijj) = the decay dynamics can be 

expressed in terms of the frequency(a;)-dependent right 
eigenvectors of W (w), the decay modes. The decay ampli¬ 
tude for each mode is determined by the overlap between 
the initial state po and the left eigenvector of W{u}) for 
the same eigenvalue, the amplitude covector. Since W 
is non-Hermitian, left and right eigenvectors are not sim¬ 
ply each others’ adjoint. However, under the very general 
assumptions stated above, we can prove that the vectors 
are still linked by the duality relation 

Wico; H, H^, M)t = -F + VW{Co-, H, (1) 

The physical consequences of Eq. (1) arise from the three 
operations it involves: (i) a constant shift by F - the 


lumped sum of constant couplings characterizing (ii) 
the fundamental transformation V, which multiplies an 
operator by the fermion-parity operator (—1)^ := 
with local fermion-number operator N; (iii) a parameter 
substitution in the original model, which constructs a 
dual model with inverted local energies H := —H, a dual 
coupling := iH'^ to reservoirs with dual chemical 
potentials ft := —p, at dual frequency ut := iT — lu*, but 
with same Hamiltonian and temperature T. The 
duality Eq. (1) therefore links in a nontrivial way the left 
eigenvectors of W to those of , which are in turn the 
right eigenvectors of W for the dual model. The physics 
behind the derivation of Eq. (1) is summarized at the 
end of the paper, and the proof is given in [47]. In the 
following, we apply the duality to the transient dynamics 
of the system depicted in Eigs. l(a,b), and show that 
it provides the key insight to interpret the remarkable 
features of the heat current. 

Transient dynamics. The system of interest is a spin- 
degenerate single-level quantum dot with Hamiltonian 
H = eN UN'^Ni, where e is the tuneable level posi¬ 
tion and U the local interaction, see Fig. 1. Here, N = 
N^ + Ni is the occupation operator, with = d^da and 
da the field operators of the dot electrons {a =t, j.). The 
dot is tunnel-coupled by a '’'k.^Cka^a + H-C. to 

a single noninteracting electrode ^ ik.aclaCka- 

Here, Cka is the field operator for reservoir electrons 
with spin cr, orbital index fc, and energy Ck^ay while 
is the tunnel amplitude. The grand-canonical reservoir 
state is p^ = , where 

= EkaCl aCka- In the wide-band limit, the tunnel 
coupling is characterized by Fo- = 2TrJ2kH^~^k)\rk,a\'^y 
assumed w-independent, with total strength F := J2a 
and F^ = F^ = r/2 in our spin-independent case. 

We consider the regime of high-temperature and weak- 
coupling to the electrode F <C T, and the response to a 
sudden shift of the level position cq —>■ e at t = 0, due to, 
e.g., a switch of the gate-voltage. Immediately after the 
level shift, the dot state equals the initial state po before 
the shift. At later times t > 0, the density operator p{t) 
of the dot obeys the Born-Markov master equation 

dt\pit)) ^W\p{t)), W := Van. W{uj), (2) 

where W = for *> J = 0; 2- Here and be¬ 

low, we write an operator x as |x) := x, and its covector 
acting on an argument • as (x|* := Tr xt* [48]. The basis 
vectors denote unit-trace physical-state operators ^ with 
0, 1 and 2 electrons: |0) := |0)(0|, |1) := \ ^a\^){'^\^ 
|2) := |2)(2|. The coefficients Wij, given in [47], repre¬ 
sent the golden-rule transition rates between the different 
electronic states (see, e.g., [34]). 


r The empty state is |0}, whereas |(t} := rf^lO) and |2} ;= d|d||0). 
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For t > 0, a standard way of solving Eq. (2) proceeds 
by expanding p{t) in the right eigenvectors of IE - the de¬ 
cay modes - and obtaining their coefficients from the cor¬ 
responding left eigenvectors - the amplitude covectors: 

|p(t))=e^‘|po) (3) 

= \z){z'\po) + e->‘b)(p'|po) + e-^=*|c)(c'|po). 

We denote the eigenvalues of W by — 7 a;, where 'jx is a 
positive decay rate for x = c,p. Since the system reaches 
a unique stationary state lim(_,.oo b(t)) = b), one eigen¬ 
value of VF is 0. 

To explicitly calculate the rates and vectors in Eq. (3) , 
one can now straightforwardly determine the left and 
right eigenvectors of the 3x3 matrix representing W. 
However, this does not give any systematic, physical in¬ 
sight into how the excitation amplitudes and rates of the 
decay modes are related. As we illustrate in the follow¬ 
ing, the duality Eq. (1) does provide such insight. In 
the Born-Markov approximation, we expand Eq. (1) to 
linear order in F, and take the limit w, zF —> iO in the fre¬ 
quency arguments. This yields ^ = —F — VTWIV, 

where T denotes the parameter substitution that con¬ 
structs the dual model: p —>■ —/i, e —?■ —e and, impor¬ 
tantly, U —>■ —U. For a known mode \x) = x with decay 
rate application of this relation to VT\x) determines 
an amplitude covector {y'\» with 

y' = {—t-)^Tx and rate = F — ly^,. (4) 

The key result Eq. (4) cross-links left and right eigenvec¬ 
tors with different rates. For our model, it allows to phys¬ 
ically relate all the rates and vectors in Eq. (3). First, 
the zero eigenvector of W is simply the trace-normalized 
stationary equilibrium state \z) = /Z of the 

quantum dot with Z{e,U,p) = Tr Its am¬ 

plitude in Eq. (3) is {z'\pq) = 1, since probability con¬ 
servation (Trp(t) = 1) for all times t requires that the 
left eigenvector associated with the zero eigenvalue is 
{z'\ = (1|. The mere existence of this left zero eigen¬ 
vector implies by duality (4) that there is a maximal [47] 
rate yp = F which depends only on the bare tunnel cou¬ 
pling, and whose mode is given by the fermion-parity 
operator: 

\p) = |(-1)^) with yp = F. (5) 

Analogously, the parity amplitude-covector is found by 
applying Eq. (4) to the stationary mode \z), giving 

ip'\ = ((-1)^a| with yp = F. (6) 


^ Since W ocF oc the dual coupling effectively 

inverts the sign of W. 


It contains Zi =Iz = z(—€, —U, —p), the stationary state 
of the dual quantum-dot model with attractive interac¬ 
tion obtained by the parameter substitution Z, i.e, by 
inverting the energies as illustrated in Fig. 1(c). There¬ 
fore, we call \zi) the inverted stationary state. Finally, the 
charge-decay rate is self-dual, yc = F—Zyc = 5 r{/+(e)-|- 
f~[e-\-U)} [33] with/=*=(e) = -I- 1)“^. Its am¬ 

plitude and mode (see [47]) are thus also connected by 
the fermion-parity and the dual model: (c'j = 2[ZZ|c)]b 

The existence of a decay mode with a rate F that sur¬ 
prisingly only depends on the bare coupling was pointed 
out in previous works ]33-37]. However, it was not un¬ 
derstood where the independence from all physical at¬ 
tributes except the “interface” property F stems from. 
The duality (4) sheds an entirely new light onto this prob¬ 
lem. It shows that the constant decay rate F is fixed by 
the fundamental requirement of probability conservation, 
via the duality based on the fermion-parity superselection 
principle. Moreover, the duality relates the complete de¬ 
cay dynamics to the stationary state \zi) of the attractive 
dual model, as seen most explicitly in (p'\ = ((—l)'^Zi|. 
It is this relation to the dual model which dictates the 
amplitudes in the heat decay. 

Heat decay. We now apply the results derived above 
to study the transient heat current out of the dot af¬ 
ter a sudden switch of the level position eg —f e. In 
the Born-Markov limit, we can evaluate this heat cur¬ 
rent [47] as Iq^t) = —dt{H — pN\p{t)), with p{t) given 
in Eq. (3), and all its ingredients obtained by the du¬ 
ality Eq. (4). This way, one not only finds the an¬ 
nounced double-exponential form of the heat current, 
Iqit) = ace~'^‘=^ OpC”'’'*’*. Importantly, the amplitudes 

ac = [e-p+^{2-Ni)U]iNo-N,)j„ (7) 

Op =U [KA’i - l)(Ao - 1) -f liPi+Po)] Ip, (8) 

can now also be expressed entirely in terms of quantities 
with a clear physical meaning: assuming that the dot is 
initially in the stationary state jpo) = b(eo)) with ini¬ 
tial level position eg, the occupation number of the dot 
No = (TVjpg) and the parity po = ((—l)^|/Og) in the ini¬ 
tial state carry all dependence on eg, as expected. How¬ 
ever, the dependence on the final level position e enters 
not only through the stationary occupation = (IVj^), 
but also through the occupation N = (fVjzi) and parity 
Pi = ((—l)'^|zi) of the inverted stationary state^ gov¬ 
erned by attractive interactions. To illustrate the con¬ 
sequence of this dependence, we plot the amplitudes Oc, 
Op in Figs. l(d,e) for a level switch eg —f e. Most promi¬ 
nently, both Oc and Op exhibit a very sharp, large change 
when tuning the final level through e — p = —U/2. As 
revealed by Eqs. (7)-(8), the reason for this jump is that 
all terms in the heat current which relate to the Coulomb 
energy dissipation, i.e., the parity amplitude Op and the 
correction cx t/ in Oc, are governed by the dual stationary 
occupation number Ni. An attractive interaction, which 
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here enters through the dual model, is well known [38] to 
force -/Vi to jump directly from 0 to 2 at e — /r = —C//2, 
avoiding Ni = 1 and keeping an even parity pi = +1. 
By contrast, a sweep of the initial level Cq — p causes 
the initial charge to traverse Nq — 0,1, 2, and the parity 
to alternate, po = 1,-1,1. This yields sharp changes of 
Up and Oc at the expected energies for a repulsive sys¬ 
tem, the two Coulomb resonances eo — /r = 0 or —U. 
We stress that while the end result for lQ{t) also fol¬ 
lows from a straightforward calculation of left- and right- 
eigenvectors of W, it does not lead to an equally concise, 
physically-motivated form of Up and Oc- Most impor¬ 
tantly, the peculiar e dependence of these amplitudes is 
only revealed by the duality. 

Another notable feature of Figs. l(d,e) is that the par¬ 
ity mode dominates Iqit) whenever it is excited: in the 
red areas in Fig. 1(e), its amplitude assumes the con¬ 
stant, maximal value Op = TU, whereas |ac| < \TU [47] 
in Fig. 1(d). This offers an interesting possibility of ex¬ 
perimentally accessing the decay dynamics of the parity 
mode. Measuring 7 p is a non-trivial task, since the parity 
mode ](—1)^) does not enter single-particle observables 
like the average charge current, which indeed decays at 
a single rate: lN{t) = -dt{N\p{t)) = {Nq - N^) 

Methods to detect 7 p have been proposed in Refs. [34, 37], 
by coupling the dot to a quantum point-contact or a sen¬ 
sor quantum-dot. Figs. l(d,e) show that the heat current 
provides a very natural and more direct way to gain this 
information. In fact, lQ{t) can be obtained using a pump- 
probe scheme by extending well-established techniques of 
mesoscopic charge detection [49, 50]. By measuring the 
rise of the stationary electrode temperature as a function 
of the delay time, one can extract the full time-resolved 
heat current. Simultaneous measurement of the charge 
current then allows to demonstrate - without fitting - 
the predicted constancy of the rate 7 p = F as well as the 
attractive-interaction signature of Op, as sketched in [47]. 

An important insight for this kind of experiments 
is provided again by the duality and is illustrated in 
Figs, l(a-c). For U T, the amplitude Op « {zi\po)UT 
essentially equals the overlap of the initial state with the 
inverted stationary state. To infer which switch excites 
the parity mode, one can thus proceed as follows: given 
the targeted final level position e, one calculates the sta¬ 
tionary state of the attractive dual system, ]zi(e)), and 
chooses the initial level position cq such that the initially 
prepared state ]po) = l•z(eo)) has the same occupation as 
|2i(e))- 

In such measurements of the heat current, one can re¬ 
solve in time that two sequentially tunneling electrons 
carry equal charge but different energy. (Only for {7 = 0, 
the heat current is simply “tightly coupled” to the par¬ 
ticle current [13], Iq (x In since Op = 0.) For exam¬ 
ple, for an initially doubly occupied dot and a switch 
expelling both electrons - the situation at the cross in 
Figs. l(d,e) - the heat current simplifies to lQ{t) « 


(e — p)lN{t) -\- {7Fe“^J’F While the excess orbital en¬ 
ergy e — /i is carried by each of the two electrons tun¬ 
neling out, the charging energy U is dissipated already 
with the first electron. Notably, when switching e to one 
of the two Coulomb resonances - e.g., at the open cir¬ 
cle in Figs. l(d,e) - the heat-current is even entirely due 
to the tunneling of the first particle, taking place on the 
shortest, temperature independent time-scale F”^. 

General duality. We conclude by discussing the main 
physical principles behind the general duality relation 
Eq. (1). The derivation of Eq. (1) is technical but 
straightforward ^ if one uses the insights established 
in [35, 36]. The main point is that in the wide-band limit, 
the best reference solution for a perturbative calculation 
of the dynamics is not the uncoupled solution — 0), 

but rather the exact solution of the coupled system in 
the limit of infinite temperature T —>■ oo. With respect 
to this solution, the propagator n(t) for the density ma¬ 
trix of the system, p[t) = Iilfi)pQ, needs to be expanded 
in terms of only part of the coupling, as a consequence 
of the fermion-parity superselection principle [35, 36]. 
What is crucial for the result Eq. (1) is that the adjoints 
of the operators occurring in this simpler expansion can 
always be expressed in parity operations [47]. This al¬ 
lows to derive order-by-order the time-propagator dual¬ 
ity n(t;i7,i7T,{^})t = e-^*VIi{t-,H,H'^,{ji})V. This 
is equivalent to Eq. (1) and holds under the general 
conditions stated in the introduction, i.e. it applies 
also to nontrivial low-temperature, strong-coupling, non¬ 
equilibrium regimes of complex multi-level fermionic sys¬ 
tems. In the commonly used expansion about the uncou¬ 
pled system = 0), the duality remains completely 
elusive. Extending the duality beyond the wide-band 
limit is challenging but seems possible using [35, 36]. 
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I. ENERGY DECAY IN A QUANTUM DOT 
A. Master equation [Eq. (2)] 

1. Golden-rule transition rates 


In the main article, the master equation [Eq. (2)] is written in abstract form to enable the fermion-parity duality 
analysis in its clearest formulation using operators rather than matrix representations. Here we give the explicit 
form in terms of occupation probabilities, which may be more familiar. To this end, one represents the density 
operator in the orthogonal basis of physical states i = 0, 1 or 2 electrons, i.e., the positive operators with unit trace 
[(l|i)]: |0) := |0)(0| (pure), |1) := \ ^cjW){'^\ (mixed), and |2) := |2)(2| (pure). Inserting a completeness relation 
XOi I0(*l/(*l*) = accounting for the non-unit scalar product (1|1) = the density operator reads \p{t)) = \i)Pi 
with Pi := (*|p(i))/(*I*) being the physical occupations of state i = 0, 1, and 2. Similarly, inserting a completeness 
relation to the left and right of W, we can express the kernel in terms of state transition rates: 

= (s-i) 

where Wij is the rate for a transition from state \j) to state |i). Substituted into the abstract equation dt\p{t)) = 
W\p{ty), we obtain the usual master equation for the probabilities 


Pa{t) 


r-2/+(e) 

/-(e) 

0 


Po{t) 

Piit) 

— -T 

2 

2/+(e) 

-f-{e)-f+{e + U) 

2/-(e + [/) 

X 

Piit) 

P2{t) 


0 

f+ie+U) 

-2f-{e + U) 


P2{t) 


The notation for the Fermi functions has been adopted from the main article: /^(e) = + 1)“^. Notably, 

the matrix elements are given by Wij ■ {j\j), and hence include the normalization of the state jj/). In the second 
column of the matrix in Eq. (S-2), this yields 1/2 times the corresponding expansion coefficients Wn. Physically, this 
factor 1/2 is caused by the spin-degree of freedom which is “integrated out” by working with the mixed basis state 
|1) = ^ point is that the rate for each elementary tunnel process between states of a fixed spin is 

|r = To-. For example, in the formula -^Pq = |r[—2/"'"(e)Po -|- f~{e)Pi], the hrst contribution reflects that the i = 0 
electron state can decay in two ways to the i = 1 electron state which comprises of two spin states; the second term 
instead expresses that a state with only z = 1 electron and a fixed spin has only one channel to decay to the z = 0 
state. 

Finally, we note that the coefficients IT 20 = IE 02 = 0, since W is evaluated only to leading, linear order in the 
tunneling strength E. The reason is that 14^20 a-nd W 20 are associated with electron-pair tunneling, and thus become 
nonzero only in the next-to-leading order (E^), see Ref. [1]. 


2. Eigenvectors of the kernel 


In the main article, we briefly sketch how to use the duality Eq. (4) in order to derive the left and right eigenvectors 
of the kernel W. However, since the main article does not provide explicit expressions for all vectors, in particular for 
the charge mode |c) and its amplitude covector {c'\, we here give a full list of all rates 'yx (minus the eigenvalues of 
W) and eigenvectors |x), (x'| for x = z,p and c. 

As stated in the main article, the kernel must conserve probability and has a unique stationary state. This implies 
the existence of only one zero eigenvalue with its corresponding left and right eigenvector: 

7.2 = 0 , ( 2 :'| = ( 1 | 

1^) = f^e) + y{e + U) + 2/+(^)r(« + um + /+(e)/+(^ + um] . (S-3) 

Note that the stationary state \z), given by the Boltzmann factor z = Tr , has been 

expressed in terms of the Fermi functions, /^, and in terms of the above defined basis states for zero, single and 
double occupation: |0),|1),|2). 
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Next, we apply the duality (4) to {z'\ and \z) in order to obtain the fermion-parity rate 7 p, its mode, and the 
amplitude covector: 

7 p = r , (p'l = [PI|z)]t = ((-l)^z,| , \p) = [(z'\XV]^ = 

1^') = ^1^) = f-(^,) + l+^, + u) + ^)|0) + 2/-(«)/+(e + m) + /■(£)/■(£ + U)\2)] . (S-4) 

The state l^:^) is the inverted stationary state of the dual model, obtained by the energy inversion operation I which 
effects H —> —H (e —>■ —e and U —> —U) and p, —> —fi. In the paper, we simply express it as the energy inverted 
Boltzmann factor Zi = Tr . Here, this energy inversion is reflected by swapping /+ •(-> f~ 

in the stationary state \z), as written in Eq. (S-3). This proves in particular that \zi) is a physical state with 
non-negative occupation probabilities. 

Finally, the charge mode |c) and its covector (c'| follow, up to a constant normalization factor, from biorthogonality 
to {z'\ and (p'| as well as \z) and |p). We find 

7 c = |[/+(e) + r(e + C/)] , (c'| = {iV|-iV.(l| , |c) = i(-l)^[|iV)-fV,|l)]. (S-5) 

Since the duality (4) holds for any exisiting mode or covector, the charge mode |c) must necessarily be related 
to the amplitude vector (c'|, and the charge rate must be self-dual. In accordance, we And 7 c = T — I"fc and 
(c'l =2 ['PTIc)]^ The additional factor 2 is the freely chosen normalization factor, which can be motivated as follows: 
by explicit calculation, one can show that {A\p) = (A|(—1)^) = 0 for any operator A that is not at least quartic 
in the dot field operators da,d]j. This means that the parity mode |p) only enters the time-dependent decay of two- 
particle quantities, implying that |c) is the only mode that affects the time evolution of single-particle observables. 
Since the charge number N is the only relevant single-particle observable of our system, the amplitude {c'\ must be 
related to the charge N, and thereby motivates the choice of a normalization factor which yields (c'| = {N\ — Af 2 (l|. 
Biorthonormality then also fixes the normalization of the mode vector |c). 


3. Why is the fermion-parity rate the highest rate? [Eq. (5)] 

In the main article, we mention that for any system that obeys the general duality relation (1), the parity rate 
given by the lumped sum of coupling constants T is the largest decay rate in the Born-Markov limit. This deserves 
some comments. We restrict our discussion*^ to master equations, i.e., problems for which only the diagonal part of 
the density operator in the energy (H) basis needs to be considered together with the corresponding blocks of the 
kernel W. In this case, the duality relation Eq. (4) dictates T to be the highest rate due to the following argument: 
for any mode x with rate there must be another mode y with 7 ^ = T — Applying the parameter substitution 
X of Eq. (4) of the main article and using that two substitutions cancel {X^ = 1), we obtain Xjy = T — Finally, 
if we assume that all rates of the system of interest as well as of its dual systemA are non-negative, the rates must 
indeed obey 0 < 7 ^, < T for all modes x. 

To argue that the dual system is expected to have non-negative rates whenever the original system has non-negative 
rates, and thus also exhibits exponential decay, let us consider the dual transformation in detail: 

• Inversion of the local energies, H ^ H = —H: For models with a finite number of orbitals considered here, the 
many-particle spectrum is finite and thus stays bounded from below after inversion of the spectrum. Thus, no 
problems with the existence of a stationary state can arise. 

• Inversion of the chemical potentials, p, —> —p: In the wide-band limit, one effectively assumes empty reservoir 
states above p and filled reservoir states below p to extend up to energy ± 00 . Inverting the sign of p does not 
change this situation. Note that we do not invert the electron energies, but only the value of p up to which the 
electron states are occupied. 

This means that the dual model comprises a finite set of discrete quantum states coupled to a flat, wide-band reservoir 
described by the same Born-Markov master equation differing only in the physical parameter values. Therefore, it 
will exhibit exponential decay whenever the original system of interest shows exponential decay, and the parity rate 
7 p = r must indeed be the highest rate. 


^ For quantum master equations in which coherences play a role, the following discussion can be extended. This is avoided here, since 
the consistent treatment of the coherences requires a careful discussion. 

^ In the Born-Markov limit, the dual system has the same coupling , cf. discussion before Eq. (4). 
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B. Heat-current 

1. Derivation of the transient charge and heat eurrent 

In this section, we derive the Born-Markov limit expressions for the time-dependent particle- and heat current. In 
particular, we explain how to obtain the excitation amplitudes of the heat current as reported in Eqs. (7),(8) of the 
main article. 

First, we employ the conservation of the total particle number N^°*‘ = to express the particle current into 

the reservoir in terms of the current out of the dot: 

lN{t) = dt (N^) {t) = -dt {N) (t) = -dt Tr [iVp(t)] = -(Ar|at|p(t)) = -{N\W\p{t)). (S-6) 

We insert the formal eigenmode expansion of \p{t)) given in Eq. (3), using the expressions for the eigenvectors given 
in Eqs. (S-3),(S-4),(S-5). In particular, we express the action of the operator N on the kernel W as 

(iV|W= (c'|W = -7c(c'| =-7c[{iV|-iV.(l|], (S-7) 

where we exploited (1|VE = (z'|IE = 0. With the trace normalization (l|po) = li we find 

lN{t) = -{N\W\pit)) = -{c'\W\p{t)) = 7 c(c'|p(t)) = 7,(iVo - (S-8) 

Here, Afo(eo) = (-^N(fo)) and Afz(e) = (Af|z(e)) are the particle numbers in the initially prepared state po and in the 
final stationary state z{e) a long time after the level shift eg ^ e (at constant interaction strength U and temperature 
T). The initially prepared state is assumed to be the stationary state with respect to the inital level po = The 

single exponential decay rate is given by the charge rate jc = 7 r[/+(e) -|- f~{e + U)]. 

The energy current into the reservoir is in general not conserved, as the tunnel coupling itself can store energy. 
Therefore, the energy flow = —dt{H)[t) out of the dot differs from the flow I^{t) = dt{H^) into the reservoir. 

Nevertheless, in the Born-Markov limit, explicit calculations show that I^{t) fs up to corrections 

oc r^, i.e., it is consistent to treat the energy current as a conserved quantity. Thus we obtain the formula mentioned 
in the main article: 

lQ{t) = dt{H^ - pN^) « -dt{H - pN) = -dt{H - pN\p{t)) = -{H - pN\W\p{t)) (S-9) 

To obtain the energy current out of the dot, we first rewrite {H\ = e(A''| + U{2\ where (A^l* = Tr N» and (2|» = 
Tr |2)(2|» = (2| • |2). Then, we use (2|p) = (2|(—1)^|2) = 1 to obtain 

{N\W={c'\W , {2\W = {2\p){p'\W+{2\c){c'\W = ip'\W+^{2-N,){c'\W, (S-10) 

where A^i(e) = (N'lzi(e)) is the average particle number in the inverted stationary state jzi) as defined in Eq. (S-4). 
With these relations, the energy current is found to be 

lE{t) = -dt{H\pit)) = -{H\W\p{t)) = - [e + 1(2 - W)H] (c'|W|p(t)) - U{p'\W\p{t)) 

= [e + i(2 - N,)U] ■ lN{t) + Ur ■ (z,|(-l)^|po)e-r‘, (S-11) 

where we have identified the charge current lN{t) = —(c'|IT|p(t)) as given in Eq. (S-8), and again used the formal 
eigenmode expansion of |p(t)) Eq. (3) as well as the eigenvectors Eqs. (S-3),(S-4),(S-5). Finally, for the heat current 
Iq (t) = IEif) — pi Nit), we obtain the double exponential form 

7Q(t) = ace-^=‘ + Ope->‘ , a, = [e - p + ^{2 - Ni{e))U] x {Ng - N,) x -f, , = UT{{-1)^Zi\po). (S-12) 

Whereas the charge amplitude Uc already has the shape reported in Eq. (7) of the main article, the parity amplitude 
here has a different form. The advantage of the form given in Eq. (S-12) is that, when analyzing the coefficients in 
the formula for \zi) in Eq. (S-4), one realizes that (—l)^|zi) fa jz,) for U ^ T. This indeed confirms our claim from 
the main article that Up fs (zj|po)17r is given by the overlap between the inverted stationary state and the initial 
state. However, what is left to show is the equivalence to Eq. (8) from the main paper, stating that 

ap/{ur) = i(W - l){No - 1) + i(pi +po), (S-13) 

with Ni, Pi depending on e, and Ng and pg depending on eg. This equivalence follows most easily when noting that the 
three operators 1, N" —1 and (—1)^ form an orthogonal basis with scalar product norms (1|1) = ((—1)^|(—1)^) = 4 
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FIG. S-1. (a) The three functions entering Eq. (S-19), namely, the energy prefactor AE, the charge rate 7c and the upper 
bound for the excess charge , as a function of the level position e with respect to the chemical potential ft for a fixed 

interaction strength U in the zero temperature limit T —> 0. (b) Upper bound for the charge amplitude, defined by 0 ®“'’ in 
Eq. (S-19), as a function of the dot level position e and the temperature T. The typical time scale is given by the coupling 
constant F. 


and {N — IjA'" — 1) = 2. Expanding the initial and inverted state in this basis, the coefficients are completely fixed 
by the average particle numbers Nq, Ni and average parities po, pi in these states: 

Ipo) = ill) + |(iVo - l)\N - 1) + ipoK-l)"^), (S-14) 

Ni) = ill) + i(iVi -1)1^ -1) + iPil(-l)^)- (S-15) 

The advantage of this basis is that the action of the parity operator - required in Eq. (S-13) - has the simple effect 
of swapping the first and last expansion coefficients, since [(—1)^]^ = 1: 

(-l)^lzO = ipill) + - l)|iV - 1) + il(-l)^). (S-16) 

Taking the scalar product with pq and accounting for the non-unit norms given above immediately shows that 
Up = 17r((— l)^Zi|po) is equivalent to Eq. (S-13). 


2. Estimate for the charge mode amplitude [Eq. (8) and Fig. 1] 

In Fig. 1(d) of the main article, we plotted the charge amplitude Uc for the specific temperature T/U = 1/10 
and pointed out that jad < . Here, we explicitly show that within the plotted window of final level positions 

—U < e — pL + < U this holds for all temperatures of immediate experimental interest, T <^U. To discuss this in 

the most efficient way, we note that according to Eq. (S-12), jad decomposes into three factors: 

|ad = |e - /r + 1(2 - Ari(e))U| x 7 ^ x |Afo - N,\ . (S-17) 

The only term that depends on the initial level position eg is the initial occupation number Nq in the modulus of the 
excess charge \Nq — Nz\. Hence, to determine an upper bound for |od as a function of only the final level e and the 
temperature T relative to the scale U, we first have to find a supremum for the excess charge \Nq — for arbitrary 
initial levels eg. Clearly, since the average number of particles on the dot is between 0 and 2, it is always true that 

0<\No-Nz\<AN^'^'i’ = ma.^{Nz{e,U,T),2-Nzie,U,T)}. (S-18) 

This means that 

^ ^ ^ AE = e-p + (S-19) 

forms an upper bound to the modulus of the charge mode amplitude jad. To get an intuition about how a®“P 
behaves as a function of the final level position, Fig. S-l(a) shows the three factors entering Eq. (S-19), that is, the 
energy prefactor AE, the charge rate 7 c, and the upper bound for the excess charge AN^'^p, in the limit T —> 0. 
Multiplying the three functions, it becomes intuitively clear that 0 < jad < < UT/2 inside the Coulomb window 

—U <e — p+^<U. However, since the limit T —> 0 is not strictly allowed to be taken in our weak coupling 
approach, Fig. S-l(b) also explicitly visualizes 0 ®“^ as a function of (e — p)/U for temperatures between T/U = 0.1 
and T/U = 0.5. Together, Figs. S-l(a,b) confirm our claim from the main paper that 0 < jad < EV[2. 
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FIG. S-2. Pump-probe scheme described in the text. 


3. Pump-probe scheme for measuring the time-dependent heat current [Eq. (S-12)] 

We briefly sketch a scheme to detect the time-dependent heat current [Eq. (S-12)] emitted from the quantum dot 
by doing a time-averaged temperature measurement in a detector electrode (D). The idea is shown in Fig. S-2 and 
involves an additional source electrode (S) to refill the dot. It requires independent voltage control over the dot 
level position and over the heights of the tunnel barriers to S and D, in order to separate the measurement from the 
refilling process. Furthermore, we assume that the electromagnetic fields arising from the voltage switches themselves 
do not lead to energy dissipation into the electrodes. Under these conditions, the detailed pump-probe sequence is 
the following: 

(a) At time t < 0, the dot is prepared in an initial state po, for example with two electrons as in Fig. S-2. The dot 
is initially isolated with respect to particle-exchange from both the source electrode (S) and detector electrode 
(D). 

(b) At t = 0, a very fast voltage switch (on the time scale l/T) brings the dot into a non-equilibrium situation with 
respect to the detector electrode by switching the quantum dot energy level cg —>■ c (by a primary gate voltage) 
such that immediately after the switch, the dot is still in its initial state po- Simultaneously, one switches on 
the tunnel coupling T to the detector D (by lowering the tunnel barrier height with an appropriate additional 
gate voltage) to enable the energy measurement. 

(c) As discussed in the main article, for f > 0 the electron tunneling will induce energy transport into the detector 
on time scales l/7c, l/7p ^ l/F, described by the heat current Eq. (S-12). 

(d) After a delay time td the measurement of the heat current is instantly shut off by isolating the dot from the 
detector (by a voltage switch that raises the tunnel barrier to D) and attaching it to the source S. 

(e) The dot energy level is switched back to its initial value e —> eg. This starts the preparation of the dot for the 
next switch. After trestart — td F the dot has again reached the stationary state pg and the system is ready 
to return to (a). 

Repeating the sequence (a)-(e) of Fig. S-2 a large number of M times, a measurable stationary temperature rise will 
result in the detector electrode due to Mx the energy that has been deposited in the time interval td- By measuring 
this temperature rise as the delay time is varied from td ^ F"^ to F"^ <C td (but respecting F"^ ^ ^restart “ ^d), one 
can extract the full time-resolved heat current lQ{t). 

It should be noted that this type of technique requires fast voltage control but allows for slow charge and heat 
measurement, and is well-established in the held of mesoscopic charge detection [2, 3]. Simultaneous measurement 
of the time-dependent charge and heat current I Nit) a-nd Iq (t) would allow a detailed experimental analysis of the 
time-resolved effects of Coulomb interaction in the heat current as predicted in the main article. 


4- Extracting the parity rate and the heat current amplitudes [Eq. (5), Eq. (7), Eq. (8)] 

In the main article, we claim that a measurement of both the transient charge and heat current after a level shift 
cg —>■ e allows to extract the parity rate 7 p and the parity amplitude Op, thereby verifying 7 p = F as well as the 
signature of attractive interactions in the amplitude in Fig. 1(e). Here, we discuss how this is possible, and also point 
out some considerations for such experiments. 
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We first outline the main idea. Using the form we obtained for the time-dependent charge current, lN{t) oc e 
and heat current, Iq {t) = a^e -|- Ope we can write 


lQ{t) = E[(x>) ■ lN(t) + ape 


(S-20) 


where E[oo) cx Oc denotes a constant with dimension energy which is to be determined experimentally. 

1. If lQ{t) and lN{t) are given, from a plot of their ratio one can extract 

E{oo) = (S-21) 

In{oo) 

Having obtained E{oo) ■ the “tight coupling” part of the heat current, one can then construct 

SlQit) ■■= lQ{t) - E{oo) ■ lN{t) = flpe-T'*’*, (S-22) 

whose value at t = 0 gives Op. In principle, repeating this for various switches eo “t e, one can, without fitting, 
completely map out the amplitude Op as shown in Fig. 1(e) of the main article. This would then demonstrate 
the signature of attractive interaction in the decay of the repulsive system. 


2. Next, plotting the normalized time-dependent quantity 



<5/q(0) 


= 7p = r. 


(S-23) 


all data points should collapse onto a single, constant value, independently of which switch eg e, which 
temperature T and which interaction strength U is considered. This would then demonstrate the exceptional 
restriction of the decay rate 7 p enforced by probability conservation, via the duality Eq. (1) based on the 
fermion-parity superselection principle. 


Our results provide an additional, independent check on the consistency of the description: since Op = UV whenever 
it is nonzero [red regions in Fig. 1(e)], dividing the values of Op by the extracted value for 7 p always gives the constant 
Coulomb energy: 


O'vhp = E- 


(S-24) 


The constancy of this ratio for all switches cq “t e well within the red regions in Fig. 1(e) is again remarkable. This 
value can be compared with the values of U determined independently by standard stationary Coulomb blockade 
linear-response spectroscopy. 

Relevant to the above discussion is the more general, time-dependent ratio with dimension of energy: 


E(t) 


iQjt) 

lN{t) 


_L °-P r-h„-yo)t 

7c {No - iV.) 7c {No - iV.) 


(S-25) 


From the expression on the right hand side, obtained by inserting Eq. (S-8) and Eq. (S-12), follow a few experimental 
considerations for the implementation of the sketched procedure: 

• Eirst, the energy levels cq and e should correspond to different charge states of the dot, such that |A^o(£o) “ 
.^z(£)| > 1, i.e. to ensure that the amplitude of In is measurable. 

• Second, one should ensure that the time scale set by the difference between the two decay rates, 7 p — 7 c, lies 

within the window of experimentally accessible times, e.g., of the pump-probe scheme suggested in Sec. IBS. 
Interestingly, our duality provides physical insight into this rate difference ^ and one finds that the timescale 
is minimal at e — /i = —U/2 where the rate-difference scales as ~ but corrections from higher- 

order tunneling processes should be expected. To avoid the degeneracy of rates for switches that extend into the 
Coulomb blockade regime - which expose the attractive-interaction features of the dual model - the temperature 
T should thus not be too low compared to U. 


^ This rate difference is always positive since jp > jc, see Sec. I A 3. In fact, this time scale is set by the charge relaxation rate of the 
attractive dual model: by our duality (4), T — 7c(e, U) = U), which equals 7 c(—e, —U) = ^ [/“(e) + /"*”(€ + t/)]. It thus makes 

physical sense that the latter effective rate is minimal at e — /i = —U{2, where it scales as ~ here, the dual model exhibits 

the electron-pair (!) charge degeneracy JVi = 0 o 2 induced by the attractive interaction, and thus allows both an empty and doubly 
occupied initial state to remain stable. 








• Importantly, the difference in amplitudes between the two terms in Eq. (S-25), required for a clear signal, 
presents no problem: as mentioned in the paper and shown in Sec. IB2, Op > 2ac, meaning that the amplitude 
of the time-dependent second term is always larger. 

Finally, we note that equation (S-25) also indicates how the charge amplitude in the heat current can be extracted: 

flc w E{oo) ■ 7c ■ (iVo - , (S-26) 

where the charge difference Nq — = /q°° dtlN{t) is, of course, available, and the required charge rate 7 c has also 

been shown to be experimentally extractable in interacting systems [4]. 


II. GENERAL DUALITY IN FERMIONIC OPEN QUANTUM SYSTEMS 

The key result of the main paper is the duality relation (1) which is straightforwardly derived in Sec. IIB using a 
formulation of open-system dynamics that was previously established in Refs. [5, 6]. However, to keep the supplement 
self-contained. Sec. IIA first introduces the required background using unified, adapted notation for the reader’s 
convenience. 


A. Previous results and notation 

The duality Eq. (1) of the main article concerns the time-evolution of the reduced density operator, p{t) = n(t)p(0), 
where p(0) is an arbitrary initial state and n(t) the propagator. As done in the vast majority of open-system 
approaches, we consider reduced dynamics starting from an initially factorized state of the total system p*°*(0) = 
p^p(O), where the density operator of the reservoirs p^ is described by the grand-canonical ensemble. Then, in general, 
the propagator is obtained by taking the trace Ttr over the fermionic reservoirs of the total system propagator, the 
solution of the Liouville - von Neumann equation = —iL*°*p*°* with = L + + L'^: 

n(t) = Tre““-^'°‘*p^ = (S-27) 

R 

The Liouvillians are the commutators with the respective Hamiltonians, = [H^, •] for x = tot, R or “empty”. Here 
- and in most expressions below - we use the same bra-ket notation [7] for Hilbert-Schmidt scalar products of the 
reservoir as used in the main article for the subsystem (quantum dot). In this case, for a superoperator S acting on 
the total system, taking the partial matrix element with respect to the reservoirs produces a new superoperator S' 
that acts only on the subsystem: 

S' = (I^IS'Ip^) := Tr5p^. (S-28) 

R 

This notation has the advantage of clearly distinguishing superoperators from the operators on which they act. 

The duality relation (1) can be expressed most conveniently by considering the Laplace-frequency representation 
p(cj) = dte®“*p(t) of the time evolution for the reduced density operator p[t) = n(t)p(0). Using that this operator 
obeys the quantum kinetic equation 

nOO 

ip{t) = -iLpit) + / dt'W{t - t')p{t'), (S-29) 

one finds 

poo ■ poo 

Pico) = n(^)p(0) , n(a;) = / dte”-*n(t) = ---, W{lo) = / dte*“‘>V(f). (S-30) 

Jo io-L-iW{Lo) Jq 

The aim of the following sections is to derive a duality relation for the propagator n(tj) that translates to the duality 
relation (1) for the kernel W(uj) by applying Eq. (S-30). 


1. Reduced dynamics in the wide-hand limit 

The expansion of the kernel W defined by Eq. (S-30) for the evolution Eq. (S-27) in powers of is the standard 
perturbation series of the Nakiyima-Zwanzig approach [8], which is equivalent [9, 10] to the basic real-time diagram¬ 
matic expansion [5, 6, 11, 12]. In the main article, we apply the duality in the very simple Born-Markov limit in which 
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W is restricted to the leading non-vanishing term in L^: 

W(t) « (S-31) 

R 

For the discussed model, this is equivalent to Fermi’s “Golden-Rule” theory. However, the new duality Eq. (1) 
applies much more generally. The starting point for our derivation is the recent observation that the expansion of W 
in the coupling Liouvillian can be drastically simplified in the wide band limit without assuming weak coupling [5, 6]. 
Based on earlier renormalization-group ideas [13], this wide-band limit formulation ^ starts from the well-known fact 
that in the high-temperature limit, the leading order perturbation theory for the exact solution Eq. (S-30) is “good”. 
More precisely, as pointed out in Refs. [5, 6], at infinite reservoir temperature, T = oo, this gives the exact solution 
in the form W(t) = Woo&{t — 0), independent of both frequency cu and the local Liouvillian L = \H, •], and results 
in strictly Markovian dynamics: 

n^(t) =^l^n(t) , W^[L^)= = (S-32) 

Here and below, we only indicate the functional dependence of the kernel W on those parts of = L + that 

are relevant for the discussion of the duality (i.e., the dependence on is not written because it is always assumed 
to be the same). The key advance made in Refs. [5, 6] was to systematically exploit this result and reorganize the 
wide-band limit perturbation theory for any finite T around this reference solution. This affords many physical and 
technical insights that remain completely obscured in the original series starting from an uncoupled system = 0). 
The most fundamental example is the duality (1) reported in the main article, which was not noted yet in Refs. [5, 6]. 
The main idea of our approach - as sketched in the main article - is that the tunnel coupling is decomposed as 
= Lf_ + Lf_ into mutually adjoint components as indicated precisely below [Eq. (S-36)]. In Refs. [5, 6], 

a renormalized perturbation series for n(a;) was derived which is exact in the wide-band limit: it is generated by 
expanding the following expression in powers of [cf. Eq. (S-61c) below]: 

U(t) = (S-33) 

In this much simplified series, one is thus allowed to drop part of the tunnel coupling, L’^ —> when simply including 
the T = oo solution iWcc{L'^) = iWrx{Lj_ + L'f) [Eq. (S-32)] into the unperturbed dynamics L ^ L + iWao- The 
derivation in Sec. HB of the key result Eq. (1) of the main article boils down to showing that this exact expression 
for n(t) obeys the duality order-by-order. 


2. Fermion-parity superselection postulate and open-system dynamics 

The mentioned splitting of the coupling Liouvillian into its parts L'^ appears naturally when using special 
quantum fields which are superoperators [5, 6]. The crucial property that makes these field superoperators such 
a powerful tool is that they explicitly incorporate the structure imposed on the open-system time evolution by the 
fermion-parity superselection postulate. We now explain at which point this postulate can be exploited. 

This superselection postulate requires [14-17] that every physical state p*°*(t) and every observable operator com¬ 
mutes at any time with the total fermion-parity operator (—1)^+^ . Thereby, quantum superpositions of pure states 
with even and odd fermion numbers = N are prohibited, as mentioned in the introduction of the main 

article. Equivalently, |p*°*(t)) must be an eigenvector with eigenvalue -|-1 of the fermion-parity superoperator, which 
applies the parity operator to both sides of its argument •: 

:= . (-1)^+^^ := [iV,.]_ , := . (S-34) 

In the Liouville-von Neumann equation = —iLf°*‘, the Liouvillian Lf°*‘ commutes with (—1)^"+^” since 

the Hamiltonian in = [H'*°‘,»] is dictated to commute with the total fermion-parity operator (—1)^+^^ to 


^ Golden-Rule rates are obtained from Eq. (S-31) when restricting it to elements connecting diagonal elements of the density operator 
(state occupations), which is sufficient for the discussed model. For a general overview of next-to-leading order corrections to the 
Golden-Rule in various formulations and the complications added by non-diagonal density matrix elements, see Ref. [10] 

^ This is equivalent to the “discrete-step renormalization” in the real-time RG framework of Ref. [12]. 

The zero in the argument t — 0 denotes an infinitesimal shift needed when Laplace transforming (5-functions. 

^ These fields are referred to as “causal superfermions” because the algebraic structure of physical expressions in terms of these fields 
expresses causality requirements. 
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preserve this requirement for all t. Hence, one finds (—• (—1)^+^ = • for any argument • appearing in any 
expression that arises in the physical time evolution of the density operator. The splitting L'^ = + LZ appearing 

in the above wide-band perturbation expansion Eq. (S-33) exploits exactly this consequence of the fermion-parity 
superselection postulate for fermionic open quantum system. Namely, we can replace the coupling Liouvillian by 

:= (S-35) 

Thus, all insights that rely essentially on Eq. (S-33) as a starting point ultimately stem from the fermion-parity 
superselection, in combination with the wide-band limit used to obtain Eq. (S-33). To our knowledge, no other 
approach makes such explicit use of this principle. Hence, many of the far-reaching consequences of fermion-parity 
superselection in open-system dynamics have not been noted so far. 


3. Decomposition of bilinear coupling 

Based on the key step Eq. (S-35) enabled by the fermion-parity superselection postulate, the coupling can be split 
as follows: 

L'^ = Ll+ LZ, LJ = Y^T 2 iGlJi[\ for q = ±. (S-36) 

21 

In Eq. (S-36), G\ and Jf are special quantum-fields ^ of the open system and reservoir, respectively They are 
superoperators, i.e., linear transformations of an operator argument - denoted by • 

Gl» =[d 2 •+q{-t)^ • (-1)’^d 2 ] , ^ l^ci •-l-g(-l)^''• (-l)^'^cij . (S-37) 

The multi-index 1 = piaiUJi labels the reservoir degrees of freedom: r^i = ± indicates particle / hole, (Ti = ± 
corresponds to spin projection j',and uji is the electron energy in the continuous wide band. To make the notation 
compact, we here use the model system of the main article, namely, a single-level quantum dot attached to a single¬ 
band reservoir as an example: 

Cr,auj = Ccrui for r/ = - and hr p = + , d^ja = d^ for rj = - and d^ for q = +. (S-38) 

This, however, does not lead to a loss of generality: one can always add additional indices to the field operators 
denoting additional orbitals, bands etc., without affecting the derivation presented here. 

The superoperators (S-37) are proper quantum fields because they obey fermionic anticommutation relations: 

[Gl Gt,]+ = S‘>’'^'S22>t, [J!, 4]+ = (S-39) 

denoting q' = —q' for g' = ±. Einally, for the particular model system, the coefficient in Eq. (S-36) reads r 2 i = 
^ 1 ( 2^1 ^ct 2 cti 1 which is independent of wi due to the assumption of a wide, flat band with energy-independent 
density of states and tunnel amplitudes in . The latter are defined as = to-i for qi = + and 

for q\ = —. The J’s express restrictions of particle and spin conservation. This admittedly dense notation 
is crucial since it captures the essence of a very large class of models to which the duality applies. The generalization 
of the assumptions and the notation is discussed below in Sec. HA 5. 


f. Expressing the adjoint in fermion-parity operations 

The duality Eq. (1) in the main article involves the adjoints with respect to the Hilbert-Schmidt scalar product 
between operators {A\B) = Tr B, i.e., (AjS'B) = (S'M|i3). Since the evolution (and arbitrary superoperators) can 


^ We emphasize that the splitting Eq. (S-36) does not correspond to a splitting of the underlying Hamiltonian into two parts, as 
one verifies explicitly by inserting the definitions Eq. (S-37) to recover (S-35). For this reason, it leads to unexpected new insights in 
open-system dynamics. 

^ See Ref. [5 and 6] for a physically motivated account of their construction. Note that in contrast to these references, we here do not 
incorporate the tunnel matrix element T 21 into the definition of Jf. 

Note that we choose the fields of the reservoir and the local subsystem to commute, [c, tf] = 0. This is always possible, as shown in 
Ref. [12]. 

This shows that the expansion (S-33) is a huge simplification relative to Eq. (S-27): since Eq. (S-33) contains only one part of the 
interaction, LAf_. in which only one type of each field appears for the system (reservoir), namely {J~)- Thus, all appearing reservoir 
fields J~ anticommute. This is not the case in Eq. (S-27), in which the additional fields do not anticommute with J~ [cf. Eq. (S-39)]. 
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be expressed [6] in the fields (S-37), we need to know their adjoints. It is readily verified from their definition (S-37) 
that ^ 


(G’)t=G| , (Jn^ = 4 (S-40) 

The adjoint thus inverts two indices. Namely, q becomes q = —q, and the multi-indices 2 = (—r/ 2 )o’ 2 j 1 = 
have inverted rj 2 ,rji relative to the multi-indices 2 = r/ 20 ’ 2 ,l = i.e., particle and hole operators are in¬ 

terchanged. An essential step in the derivation of the duality relation (1) is to express the index inversion of the 
adjoint (G|)^ in Eq. (S-40) using parity operations instead. This can be done by composing three relations established 
in Ref. [5], and works as follows: 

- First, consider the operator-adjoint transformation JC, which maps an arbitrary operator • into its adjoint 

/C|.) = |.t). (S-41) 

This operation is of fundamental importance, since K,p = p must hold for any valid density operator [cf. 
Eq. (S-47)]. By definition, it inverts the rj index in the multi-index 2 of the ordinary fields: lCd 2 = d 2 - The 
transformation of the fields G 2 with JC essentially does the same: 

/C(G«)/C = gG?(-l)^'^. (S-42) 

This follows directly from the definition (S-37), the essence of the derivation being that due to JC{AB) = , 

the operator ordering is reversed in Eq. (S-37), effectively swapping its two terms resulting in a sign g = ± and 
a parity operation 


(-1)^". = (-1)^.(-1)^. (S-43) 

Below, we will use that (—1)^” anticommutes with all fields G 2 and thus counts their fermion-parity (see note 
with Eq. (S-44)) and squares to one, [(—1)^^]^ = 1. 

- Now consider the other, independent parity transformation 

V, ■= (-1)^. (S-44) 

featuring in the central result Eq. (1) of the main article and not to be confused with (—1)^” = (“1)'^ • (“1)^ 
appearing in Eq. (S-42). The transformation of the fields with V only inverts the <? = ± index: 

V{Gl)V = -Gl- (S-45) 

This also follows directly from the definition (S-37): essentially, by left-multiplying with (—1)^ and anticom¬ 
muting it through di», the sign of the first term in Eq. (S-37) is inverted, which is equivalent to inverting q and 
an overall sign. 

Composing Eq. (S-40)-(S-45), we see that physical expressions involving the adjoint - which are always of the type 
[/CLT/C]t (see below) - can always be expressed in terms of linear fermion-parity operations: 

[ICGllC]'’= qV{-l)^"'GlV, (-1)^" = (-1)^. (-1)^, p. = (-l)^.. (S-46) 

Equation (S-46) is crucial in proving the duality relation in Sec. IIB. 


^ See Eq. (S-80) in App. IID land also Eq. (57) from Ref. [5]. 

^ One should note that the superscipt index q- = ± in Eq. (S-37) is not related to physical real particle / holes: this is rather regulated 
by the separate index 77 = ±. 

See Eq. (S-81) in App. IIDl below and also Ref. [5] [Eq. {G- 6 )]. 

See Eq. (S-82) in App. IIDl below and also Ref. [5] [Eq. {E-7)-(E-8)]. 
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5. Applicability to fermionic open-system models 

In the following, we restrict our attention to the simple model of the main article as mentioned with Eq. (S-38). 
However, the notation is so general that without modification, the derivation of Eq. (1) goes through under much 
more general conditions: ® 

• System: Nowhere did we use special properties of H, except that it preserves fermion-parity, as it should by 
the postulates of quantum mechanics. It describes the local evolution of a finite set of discrete many-fermion 
quantum states due to any finite number of orbitals. The Hamiltonian H can be of arbitrary complexity, including 
all possible multi-fermion interactions (not just “charging energy” as in the Anderson model). However, we stress 
that we restrict our attention to purely fermionic models H, i.e., it does not contain local bosonic degrees of 
freedom (e.g., as in the Anderson-Holstein model). 

• Reservoirs: There can be multiple reservoirs, each with multiple hands. Extending the multi-index 1 = lyiaioji 
by additional band or reservoir indices takes care of this. 

• Coupling: We do not have to specify any property of the coupling , except that it is bilinear and energy- 
independent. The special simplification possible in the heat-current problem analysed in the main article, namely, 
that it conserves spin-projection and total spin length, are thus irrelevant. Each band in each reservoir can have 
completely general and independent couplings to each spin-orbital of the system (within the bilinear, flat, wide¬ 
band limit). The only change is that in Eq. (S-54), the last expression on the right for E will be modified to a 
different constant by summing over more bare coupling constants. However, it remains an “interface”- constant. 


B. Derivation of the fermion-parity duality relation [Eq. (1)] 

We now write out the derivation of the duality relation Eq. (1) in the main article. As mentioned there, it follows as 
a straightforward consequence of the prior results established in Refs. [5, 6] and reviewed above, even though it requires 
some notation. The hrst observation to make is that this duality features the adjoint operation f. This antilinear 
operation, when applied to physical superoperators, is not tied to any conservation of fundamental properties: for 
an open system, the dynamics are not unitary, i.e., p{t) ^ U{t)p{0)U{tf [for some unitary U{t)], which implies 
that n(t)f n(— t). An important ingredient of the derivation is to express the effect of taking the adjoint of the 
propagator in terms of operations that do relate to symmetries that are preserved by open system time evolution. 
One such fundamental symmetry is the selfadjointness of the reduced density operator p{tf = p(t) that must hold 
for all times t. With p{t) = n(t)p(0) and the operation K of taking the operator-adjoint defined in Eq. (S-41), this 
requirement is 

/Cn(f)/C = U{t). (S-47) 

This is, of course, satisfied by both the exact general propagator [Eq. (S-27)], but also by the exact wide-band limit 
propagator (S-33) [and (S-49) below] required here. The derivation of the duality relation Eq. (1) exploits this 
symmetry by studying the combined linear operation [K. • Afjt on physical superoperators instead of the antilinear 
adjoint f. Just as the adjoint, this operation reverses the order of products: 

[KBAK]^ = [/CA/C]t[/CB/C]t. (S-48) 

By applying [K. • Kf to the exact wide-band evolution propagator (S-33), 

n(t) = (S-49) 

we can derive a duality relation for n(t) [see Eq. (S-76)] that, by virtue of the solution to the generalized quantum 
master equation, Eq. (S-30), is equivalent to the duality relation (1) shown in the main article for the kernel W in 
the Laplace frequency representation [Sec. HB3]. The benefit of this procedure is that while the duality is most 


^ We note that U(l) symmetry, i.e., particle number conservation is allowed to be broken in the system and by the coupling, but we used 
particle conservation in the reservoirs [Eq. (S-52)]: ^21 ^ treat mean-field superconducting electrodes, one may consider 

lifting this constraint, but the energy-dependent density of states in this case does not allow Eq. (S-33) to be derived from Eq. (S-27). 
This is an interesting future topic. 

^ One readily verifies to be satisfied by the general propagator Eq. (S-27): using fCiL^fC = iL^ for fc = R, T or ”empty“, one finds 

— Trpi + — II(t). In Ref. [5] [Appendix G], the expansion (S-61c) of the exact wide-band limit propagator 

Eq. (S-33) was explicitly shown to satisfy it order-by-order. 
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straightforwardly proven in the time-representation, the applications of the duality become more obvious in Laplace 
frequency space, e.g., if one considers the zero frequency Born-Markov limit as in the main article. 

In the following, we proceed in two physical stages in the spirit of Refs. [5, 6]: first, we consider the simple inhnite 
temperature T = oo limit [Sec. IIB 1]. Then, we use this to simplify the exact time-evolution for finite temperatures 
T < oo [Sec. IIB 2]. 


1. Infinite temperature T = oo limit: exact reference solution 

We now hrst consider the simple solution at T = oo around which the expansion of (S-49) is organized; 

HocW = L^=L + iW^. (S-50) 

Since Boo is a physical limit of a valid time-evolution, it must preserve hermiticity. Indeed, it follows from the 
definition of the kernel Woo in Eq. (S-32) and L* = [H, •] that 

KL1C = -L , ICWoolC^Woo ICiLooKi^iL^ , K:n^{t)K: = (S-51) 

Below, we use that this also holds for the inverse noo(t)“^ = noo(— t) since Eq. (S-51) holds for all t, including t < 0. 

We first derive the duality (1) for this T = oo case. To this end, we express ^ VLoo, given by Eq. (S-32), in the fields 
Eq. (S-37): 

Wo. = -J2i^2’i'G+G^,, :=27Tj2T2'2Tri62i. (S-52) 

2'!' 21 

Here and below, the unprimed (primed) multi-indices 1, 2 (1', 2') refer to the reservoir (local subsystem). 

At first sight, no non-trivial duality relation seems to arise in this limit, as the kernel is simply self-adjoint: 
Wjo = Woa- However, realizing that selfadjointness does not extend to hnite temperatures T < oo, ^ we instead 
exploit the above discussed symmetry relation = [AHTooAIJf to express the adjoint in parity operations via 
Eq. (S-42), and then find a non-trivial duality that does have a finite-temperature extension: by inserting Eq. (S-46) 
into Wjo = [/CWoo/C]^ and using Eq. (S-43), the fields in the expansion (S-52) of Woo are permuted [Eq. (S-48)]. To 
restore their original order, we use the fermionic anticommutation relation, yielding 

WooiL'^y = -^p[Y,^2'vi-l){-t)^"G^,{-t)^"G+]v = -^v[Y,^2'vG^,G+]p 

2'V 2'1' 

= (<52-1' - G+G^,) P = -r- PWoo{L^)P. (S-53) 

2'V 

This forces a shift by the real positive number 

r = ^|r2'i'52'i'(S-54) 

2'1' V 

simplifying to T = To- for the simple model of the main article, and an overall sign change. Using that Woo is 
quadratic in according to Eq. (S-32), we can absorb this sign into a dual coupling . With this, we obtain 

the duality (1) of the main article for the limit T = oo: 

WooiL'^y = -r + PWoo{L^)P. (S-55) 

Its main usefulness for the following lies in the duality that it implies for noo(t): its adjoint can be related to a dual 
T = oo evolution noo(t): 

Uoo{ty = [/Cnoo(t)/C]t = ^ e-^*Pnoo(t)L’. (S-56) 


^ See Eq. (38) and (65) of Ref. [6]. 

^ The relation = Woo relies on the fact that both types of fields, and G~, appear in its expression (S-52). Taking the adjoint, 
one permutes two fermion fields and but also interconverts them, (G^,)^ = Gj. Using the general property T 21 = T 21 , the original 
form is then restored. As a result, left and right eigenvectors of Woo are mutually adjoint, as for Hermitian Hamiltonians of closed 
systems. However, for finite temperature T < 00 , no corresponding relation holds because the structure of the finite T corrections is 
very different: they contain only G^ fields, as we will see below [Eq. (S-61c)]. 
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The dual evolution ncxD(t) is generated by a dual frequency-independent dissipative Liouvillian Loo'- 

noc(t) = noo(t;i,T'^)=e-*^“*, Loo = L-iWooiL'^) = L + iWoo{L'^)- (S-57) 

Here and below, dual quantities are indicated by an overbar. ^ They determine the dynamics of an underlying dual 
model with the same reservoir Liouvillian (L^ = is never written), but with inverted energy signs in the local 
Liouvillian L and couplings multiplied by i, recalling a Wick rotation: 

L ■= -L, ■-= iL^. (S-58) 

What remains to be proven in the following is that the duality (S-56) with the dual Liouvillians given in Eq. (S-57) 
extends to finite temperatures without making any further approximations. 


2. Finite temperature T < oo: exact wide-band evolution 

The duality for the finite T wide-band evolution (S-33) follows most clearly when working in an interaction-picture 
which is defined in the usual spirit by “extracting” the simple reference dynamics from the exact evolution: 

n(t) := noo(t)n(t,o)/. (S- 59 ) 

Notably, in our case, this extracted evolution is already dissipative, as the reference solution (S-50) applies to the 
coupled system in a certain physical limit. In Eq. (S-59), the time-argument (t, 0) [instead of just (t)] serves to indicate 
that the interaction picture is defined with respect to the initial time 0 which will matter below. Since Hoo (t) preserves 
self-adjointness for all t [Eq. (S-51)), so does n(t,0)/ = KIi{t,Q)iK . The derivation of the duality therefore reduces 
to first expressing n(t, 0)| = [Afn(t, 0)/Af[t in a dual evolution analogous to the T = oo case [Eq. (S-56)), and then 
transforming it back to the Schrodinger picture using 

[/cn(t)/c)+ = [/cn(t,o)7/c)t [/cnoo(t)/c)^ = e-^* [/cn(t,o)7/c)t vnoo{t)r. (s- 60 ) 

Here, we used Eq. (S-48) and the T = oo duality (S-56). 

In this interaction picture, the perturbation expansion for the wide-band evolution n(t, 0)7 takes the familiar form: 
n(t, 0 )7 = J'o jp^) (S-61a) 

OO „ 

= '^{-if I dtk-.-dti (ll^|L(]((4)...L[^(ti)|p^) (S-61b) 

OO p 

= E E / dtk-.-dh n>k---TviG+{tk)...G+{h){l^\J-{tk)...J^{h)\p^), (S-61c) 

k even 

where T denotes the time-ordering. In the last line, ^ Eq. (S-61c), the order k is limited to even integers by the 
fermion-parity superselection In Eq. (S-61c), the Heisenberg field superoperators with respect to t = 0 are defined 
in the standard way for the reservoirs, whereas for the subsystem, we use the exact T = 00 dissipative dynamics 
(S-50): 


Jf(t) Gf,(t) := noc(-t) noo(t). (S-62) 

The time-dependent perturbation expansion (S-61c) was established in Ref. [6), based on earlier Laplace-space for¬ 
mulations [5, 12j. This perturbation expansion is specific to the wide-band limit: it relies on the assumption that the 
couplings Tpi do not depend on the reservoir frequencies u)i contained in the multi-index reservoir sums in Eq. (S-61c). 


^ The notation used in this supplement differs from that in Refs. [5, 6] where the overbar has a different meaning. 

^ To keep the notation simple, k as an index on T^/k or G^, denotes the k-th multi-index, whereas is the k-th time argument. Similarly, 
the first sum is over all possible even values of k, while the second sum is a sum over multi-indices. 

It is not Wick’s theorem which requires this: at this point we do not yet assume noninteracting electrodes. 

^ Here, the time argument suffices to indicate the interaction picture. 
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As discussed in Sec. IIA 5, each such multi-index sum is a sum over all indices necessary to characterize the system. 
For the model considered in the main article, it explicitly reads 

OO 

Since the fields of the local subsystem by definition do not depend on reservoir frequencies a;;, the reservoir 

frequency dependence in the wide-band limit is entirely contained in the reservoir fields. Exploiting this simple fact 
and furthermore using that only one type of reservoir field (Jfi) occurs in (S-61c), the finite temperature duality in 
the dissipative interaction picture is found to take the elegant form 

[KTl{t,Q)iK:]^ =VYV{<d,-t)iV , n{0,-t)i :=n{0,-t-L,L'^,fi)i. (S-64) 

In addition to the dual quantities Eq. (S-58), this involves the dual chemical potentials: 

ft := —/i (S-65) 

We derive Eq. (S-64) from two basic observations: 

1 . The adjoint f combined with operator-adjoint 1C* 1C can be expressed entirely in terms of parity operations. The 
established result Eq. (S-46) shows how this works for Schrodinger picture fields. Using the T = oo duality 
(S-56), we can extend this elegantly to interaction-picture fields (S-62): 

[/CG«,(t)/C]t = Glit)-.= n^{-t)G'l,n^(t), (S-66) 

where on the right hand side, we now have G^,, the field G\, in the interaction picture with respect to the dual 
T = oo evolution [Eq. (S-56)] generated by Loo [Eq. (S-57)[, as indicated by the overbar. Applying [1C • /Cjt to a 
sequence of an even number k of fields of the same type g = -h, their order is reversed [Eq. (S-48)[, their times 
are inverted, and a total sign (—1)^/^ appears: 

[/CG+(4)...G+(G)/C]t = (-i)^~G+(-G)...(-1)^”G+(-4) = (-1)'=/2G+(-G) ... G+ (- 4 ). (S-67) 

The sign arises by pairing up all parity superoperators (—1)^" and using [(—1)^"]^ = +1 [Eq. (S-46) ff.]. 

2. Under the reservoir frequency integrals and reservoir average, we are free to simultaneously invert the ordering 
of the fields J~, their time-arguments and all chemical potentials. We can thus write 

^L’fc„.i(lR|J-(4)... (1^1 Jr(-4)... , (S-68) 

k..l k..l 

where Fk...i denotes any expression depending on all indices contained in the sequence of multi-indices k.. .1 
exeept for the frequencies ujk,...aii. It requires two considerations to see this. First, inverting the order¬ 
ing of the fields gives a sign, since fields of the same type q = — anticommute [Eq. (S-39)], and the sign 
equals — (— 1)^/2 fQ,. gyg^ Second, by inverting all times, we transform Jf(t) = 

to Jf{—t) = J®. The sign in the exponential function can be reverted by changing the sign of all 

frequency integration variables, which is possible since we integrate from —00 to 00 . However, the inverted 
frequency signs then also appear in the Wick expansion ® of the correlator (l^j J^T ... given by the 

sum of all permutations of pair contractions (“1)^ n(n () since each Wick contraction 

depends on the frequency through ^ 

(li'l Jr|p«) = ^21 tanh . (S-69) 


^ Note that in perturbation theory renormalized by continuous RG transformations, the vertex operators G do pick up a frequency 
dependence. For this reason, it is an interesting open question if our derivation can be modified to show that the RG preserves our 
duality or breaks it. 

^ Note that non-trivial dependence on all other reservoir quantum numbers is allowed for. 

Notably, the sign ( — 1)^ appears three times, and thus survives. This sign leads to the imaginary factor in the duality . The 

sign appears once in Eq. (S-67) and twice in the derivation of Eq. (S-68) where it cancels. Up to Eq. (S-67), all expressions are valid 
for arbitrary initial nonequilibrium reservoir state (i.e. possibly interacting). It is only in the last step establishing Eq. (S-68) that 
we used that an equilibrium reservoir state is determined by the noninteracting reservoir Hamiltonian This crucial property of 
the noninteracting, equilibrium reservoirs ultimately allows the adjoint to be expressed in forward time-evolution of some dual model 
in the wide band limit. 

Schematically, (-1)^” G = (-l)GG = -GG. 

^ See Eq. (60) and Appendix A of Ref. [5]. 

^ See Eq. (S-83) in App. IID 2 below and Ref. [6] [Eq. (101)]. 
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Hence, to leave all pair contractions unaltered apart from an overall sign, the frequency inversion must be 
accompanied by the inversion of all chemical potentials fi. Altogether, this yields 

(S-70) 

For fc/2 pair contractions (even k), one obtains a sign (—1)*/^, which however cancels the sign due to the 
reordering of the reservoir fields in Eq. (S-68), and thereby verifies the latter equation. 

The duality (S-64) follows immediately from the general structure of the wide-band limit perturbation expansion 
(S-61c) and the two observations pointed out above: taking the adjoint using Eq. (S-67) gives Eq. (S-71a) and 
corresponds - apart from a parity transformation V - to inverting all times as well as the order of the operators in 
which they appear. Using Eq. (S-68) to obtain Eq. (S-71b), we recover a forward time-evolution, starting at —t and 
ending at time 0: 


[/Cn(t,0)7/C] 


CXD „ 

= !!(-*)'' H / dtk...dti Tk'k...Tv^V{-lf/^G+{-ti)...G+{-tk)V{t^\Jk{tk)...J^{ti)\p'^) (S-71a) 

fc=0 k'k...Vl-^ t>tk>...>ti>0 

k even 

oo ^ 

= / dt^...dtk Tui...rfc,fci'=G+(-ti)...G+(-4)(l«|Ji-(-ti)...4-(-4)|p^)P 


k=0 
k even 


oo « 

= dh...dtk {t^\Ll{-t,)...Ll{-tk)\p^)V 

= = VU{Q,-t)iV 


(S-71b) 

(S-71c) 

(S-71d) 


Note that the dual time evolution H is defined with respect to the same reservoir Liouvillian L^, but a dual reservoir 
density operator p^ = in Eq. (S-71b) with inverted chemical potentials p. 

Clearly, n(0, —t)/ and n(t,0)/ are related by a time-translation (not: time-reversal!), since our total system is 
time-translation invariant. In the Schrddinger picture, the time evolution is equal in both cases, but in the interaction 
picture, a time-shift results in a shift of the reference time: 

J){-ti) = J({t - U) (S-72) 

Gi{-ti) = n^{ti)Gin^{-ti) =n^{t)Gi{t-ti)n^{-t). (s-73) 

Since the transformation Eq. (S-72) drops out under the reservoir average in the step (S-71b), we obtain: 


[/cn(t, o)7/c]t ^(o, -t)iv = vn^{t)n{t, Q)in^{-t)v. (s-74) 

The Schrddinger picture dual evolution is defined, analogously to Eq. (S-59), by extracting the dual trivial T = oo 
evolution: 


n(t) = noo(t)n(t, 0)/. (S- 75 ) 

Finally, inserting this and Eq. (S-74) into Eq. (S-60), and using preservation of hermiticity, H = JCIiK. [Eq. (S-47)], 
as well as = 1, we arrive at 

n(t)t = [/cn(t)x:]t = e-^*vtt{t)v. (s-76) 

This is the general fermion-parity duality for the propagator in the time-representation. 


3. Kernel duality infrequency representation [Eq. (1)] 

Equation (S-76) essentially completes the derivation. It remains to convert it [cf. Eq. (S-30)] into the equivalent 
kernel duality Eq. (1) reported in the main article, which is more convenient for applications. With our definition of 
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the Laplace transform, /(w) = /q°° f {t), applying dte^'^ * to the duality (S-76) yields the propagator duality 
in the Laplace-frequency representation: 


rro' r rT Mt (S-30) —l 

[n(aj;L,L ,/r)]T = -^^ 

= -^ = Pn(a;;L,L ,^)P. (S-77) 

The second line is the transform of 'Pn(f)'P at shifted complex frequency ui = —ui* + iV. Using Eq. (S-30) and L = —L, 
we can thus directly read off the duality between the kernels, the central result Eq. (1) of the main article: 

[lU(w; L, L^, At)] ^ = -T + VW{oj: L, L^, fj,)V (S-78) 

In the main article, the functional dependence on the Hamiltonians was written. Here, this dependence enters through 
L = [.H,»], and 


C. Duality as a consistency check on calculations [Eq. (1)] 

As mentioned in the discussion of the duality (1) in the main article, the right hand side of Eq. (S-78) can always 
be calculated if the kernel W has been obtained, either exactly or approximately, as a function of the parameters 
specifying the model, i.e., the parameters appearing in H, , and the p’s in p^. One simply substitutes these 
parameters into W (iw, H, , p) such that H —> —H and —> iH'^, respectively, and inverts all p’s in this known 

function to obtain ,p). For example, for the Anderson model in the wide-band limit with spin-dependent 

tunneling, E = the kernel duality reads: ® 

[W{UJ-e, U, r„ p)] ^ = -T + VW{ir - w*; -e, -U, -E,, -p)P. (S-79) 

Since we have proven the duality in a general context, it can thus be exploited for any model that satisfies our 
assumptions [see Sec. HAS] when performing some perturbative or nonperturbative calculation of the kernel. 

Conversely, it is easy to check explicitly if a calculation actually preserves the exact relation Eq. (1). For example, 
we have verihed that in the limit of vanishing interactions but for arbitrary coupling E relative to T, the duality 
is satisfied by the exact, wide-band limit solution for the nonequilibrium, non-Markovian dynamics of an arbitrary 
multilevel Anderson model. The exact reduced density operator for this problem, i.e., including all multiparticle 
correlations on the subsystem (which factorize in this limit), was computed in Ref. [6]. The result for the exact 
propagator n(t), written in the time representation, was written in a compact form that most explicitly satishes the 
duality Eq. (S-76) when the above mentioned substitutions are performed. 


D. Collection of some relevant derivations 

To keep this supplement self-contained, we finally write out and comment on some of the algebra leading to the key 
formula (S-46) for converting adjoints into parity operations, and to the reservoir correlator (S-69). We emphasize 
that these results were already established in Ref. [5] and are rederived here for the reader’s convenience. However, 
since they provide the crucial starting point for the novel derivation in Sec. H B of Eq. (1) of the main article, we here 
calculate them in a form that relies less on the context of Ref. [5]. 


1. Superoperator adjoint [Eq. (S-jO)], operator-adjoint [Eq. (S-42)[, and parity [Eq. (S-45)[ 

Relation (S-46) was obtained by simply composing (S-40),(S-42) and (S-45). The latter three relations all follow 
directly from the definition of the field superoperators (S-37). 


Note that the couplings also feature in the frequency argument through their lumped sum F = Fo-: here, the Fcr’s do not get inverted 
(they are not model parameters here, but part of the imaginary frequency shift). 
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a. Proof of Eq. (S-fO) (adjoint): For arbitrary operators A and B, the adjoint (G^)^ t)® written as 

{A\{Gl)^B) {GIA\B) ^Tr [{d^A + A{-1)^B] (S-80) 

= j^Tr[{AUl + qdU-l)^AH-l)^}By''="'>^Tr[A^d^B-q{-l)^B{-l)^d-^}] = (^|G?i3). 

We used (1) the dehnition of the adjoint of superoperator S, (AjS'fi?) := (5^15); (2) the definition of the field 
superoperator, Eq. (S-37); (3) the definition of the overbarred multi-index 2, = d^ [Eq. (S-41) ff.]; (4) the cyclic 

trace property and, in the second term, the property d 2 (—1)^ = —(—l)^d 2 of the fermion-parity operator. The 
proof for (Jf)t is equivalent. 

b. Proof of Eq. (S-42) (operator-adjoint): Consider the action on an arbitary operator ji?) = B: 

JC {GKJCB)) /C(GlSt) j=JC [d2B^ + qi-t)^B^-lfd2] ® ^ [Bda + qd2(-tf B{-1)^] 

- 72 [{-'^)^B'{-l)^d2 -b qd2B'] qGlB' ® gG?(-l)^”s. (S-81) 

Here, we employed: (1) the definition of /C, KB := i?l; (2) the definition of the field superoperator Eq. (S-37); (3) the 
definition of K and K{AB) = B^A^ = {KB){KA). (4) We temporarily set B' = (—1)^ = (—and thus 
B = (—1)^, (5) extract the sign q using g = ±1 to recover the definition Eq. (S-37), (6) and finally restore 
B again. 

c. Proof of Eq. (S-jS) (parity transform): Consider the action on jB) = B, where B is an arbitary operator: 

V{Gl)VB (-1)^{g«{(-1)^b}} ® ■^(-l)^{d2(-l)^B-bgB(-]l)^d2} 

I ^ [-daB -b g(-l)^B(-l)^d 2 ] = -Gf with q := -q. (S-82) 

We applied (1) the definition of P as left-multiplication by the parity operator, P» = (—1)^», (2) the property 
[(— 1)^]^ = 1, (3) and the property d 2 (—1)^ = —(—l)^d 2 of the fermion-parity. 


2. Reservoir correlator (S-69) 

Finally, we write out the calculation of the single (!) reservoir correlation function (S-69) that appears in expansion 
(S-61c) of the wide-band limit propagator (S-33). Although its expression is well-known we write out the calculation 
to show how it arises from the definition of the fields (S-37). Inserting this definition, one calculates Eq. (S-69) as 
follows: 

TrJ^-Jfp^ Tr^J 2 -{[ci,p«]+} ® iTr[c 2 , [ci,pi']+]_ W[c 2 ,ci]_pi' 

- '^SlTr [Sr,2,+ [cUiV2’(^^2UJt,]- + - ijp'*’ 

- ^ 2 iV 2 (2/((w2 - p)/T) - = S 21 tanh (^ ^ 2 (^ 2 ^ ^ 2 ) ^ ^ 

Here, we used the following simple facts: (1) the definition of fields (S-37) and that the reservoir state - by the 
superselection postulate - must have even electron parity, (—1)^ p^(—1)^ = p^ , i.e., it has no unpaired electron 
operators; (2) one “excitation” creates odd-parity operators, (—1)^ [ci, p^]+(—1)^ = —[ci,p^]+ using (—1)^ ci = 
—ci(—1)^ ; (3) the trace is cyclic; (4) the expectation is nonzero only if 1 = 2 = (—P 2 )cr 2 W 2 and 2 = 1)202022 with 
P 2 = ±; (5) commutators are antisymmetric; (6) Ttc]^^^^Ccy 2 LJ 2 P^ — fi^) = (e^ “ 1)”^ with x = (w 2 “ P 2 )/T\ (7) 

R, 

2f{x) — 1 = tanh(x/2) and 772 tanh(x/2) = tanh(p 2 a^/ 2 ) for 772 = ±- 
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